Available  online  at  www.sciencedirect.com 


ELSEVIER 


*V  ScienceDirect 


Journal  of  Power  Sources  162  (2006)  444-459 


JOURNAL  Of 


www.elsevier.com /locate /jpowsour 


Effect  of  flow  and  temperature  distribution  on  the 
performance  of  a  PEM  fuel  cell  stack 

Jaewan  Park,  Xianguo  Li* 

Department  of  Mechanical  Engineering,  University  of  Waterloo,  Waterloo,  Ontario,  Canada  N2L  3  G1 

Received  16  April  2006;  received  in  revised  form  14  July  2006;  accepted  18  July  2006 

Available  online  1  September  2006 


Abstract 

A  non-isothermal  stack  model  has  been  developed  to  analyze  the  effects  of  flow  variance  and  temperature  distribution  on  the  performance  of  a 
polymer  electrolyte  membrane  (PEM)  fuel  cell  stack.  The  stack  model  consists  of  the  flow  network  solver  for  pressure  and  mass  flow  distributions 
for  the  reactant  gas  streams  and  cooling  water,  and  the  heat  transfer  solver  for  temperature  distribution  among  the  cells  in  the  stack,  as  well  as  the 
fuel  cell  model  for  individual  cell  performance.  Temperature,  pressure  and  concentrations  of  fuel  and  oxidant  are  the  most  important  conditions 
for  the  fuel  cell  operation.  In  this  work,  pressure,  temperature  and  concentration  distributions  are  determined  incorporating  the  individual  cell 
performance  with  the  minor  losses  in  stack  flow  network  accounted  for.  The  results  indicate  that  the  effect  of  temperature  is  dominant  on  the  cell 
voltage  variance  when  the  flow  variance  is  small  for  sufficiently  uniform  distribution  of  reactant  flow  among  the  cells  in  the  stack.  Sufficient  flow 
uniformity  can  be  achieved  by  a  large  manifold  that  reduces  the  cell  active  area,  or  a  small  flow  channel  diameter,  the  latter  may  result  in  excessive 
pumping  power  for  the  anode  and  cathode  gas  streams.  The  manifold  and  flow  channel  diameters  were  optimized  considering  stack  performance 
and  reactant  stream  pumping  power  requirement.  It  is  further  shown  that  the  flow  and  temperature  distribution  have  a  different  influence  on  the 
stack  performance,  and  a  judicial  matching  of  their  distribution  can  provide  the  ideal  uniform  cell  voltage  distribution.  An  optimal  combination  of 
the  flow  and  temperature  distribution  along  the  stack  yields  the  optimal  stack  performance. 

©  2006  Elsevier  B.V.  All  rights  reserved. 
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1.  Introduction 

The  greenhouse  effect  is  currently  considered  as  one  of  the 
most  pressing  environmental  problems  caused  by  the  use  of  fos¬ 
sil  fuels.  Such  environmental  concern  brought  about  the  strong 
demand  for  clean  power  generation  and  fuel  cell  is  receiving 
increasing  attention  as  a  clean  power  source,  at  least  at  the  point 
of  use.  With  its  competitive  power  density  and  high  efficiency  as 
well  as  zero  emission,  the  polymer  electrolyte  membrane  (PEM) 
fuel  cell  has  emerged  as  a  promising  alternative  to  internal  com¬ 
bustion  engine  [1].  Cells  are  connected  in  series  to  form  a  fuel 
cell  stack  in  order  to  satisfy  the  power  demand  of  the  practi¬ 
cal  applications.  The  performance  degradation  is  known  to  be 
caused  by  the  unequal  supply  of  fuel  and  oxidant  which  results 
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in  variations  of  each  cell  flow  conditions  for  isothermal  tempera¬ 
ture  distribution  throughout  the  stack  [2].  In  reality,  temperature 
variation  along  the  stack  has  a  significant  impact  on  stack  perfor¬ 
mance,  which  is  known  as  the  thermal  management  that  remains 
one  of  the  critical  issues  for  PEM  fuel  cells.  Therefore  the  stack 
performance  must  be  analyzed  based  on  the  actual  conditions 
present  for  each  cell  in  the  stack,  which  are  influenced  by  the 
design  and  operating  conditions  of  the  stack.  The  optimization  of 
stack  design  and  operating  parameters  is  necessary  for  improved 
performance  and  longer  lifetime. 

A  considerable  body  of  literature  exists  for  the  modeling  and 
simulation  of  a  single  PEM  fuel  cell,  as  summarized  in  a  recent 
review  article  [3].  Most  of  these  analyses  and  modeling  efforts 
can  be  categorized  into  three  categories:  the  first  is  the  empirical 
correlations  for  the  cell  performance  with  limited  applicability 
and  lack  of  predictive  capability  [4] ;  the  second  is  the  detailed 
CFD-based  simulation  [5-23].  These  multi-dimensional  CFD- 
based  models  impose  significant  computational  penalty  so  that 
the  models  are  limited  to  either  over-simplified  cell  geome- 
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Nomenclature 

Aceii  active  cell  area  (m2) 

Af  flow  path  area  (m2) 

Cf  wall  friction  coefficient 

Ci  concentration  of  species  i  (kmol  m-3) 

Ci  laminar  wall  friction  coefficient 

Cmin  minimum  heat  capacity  rate 
Cp  heat  capacity  at  constant  pressure  (J  K-1) 
dh  flow  channel  hydraulic  diameter  (mm) 

Dh  manifold  hydraulic  diameter  (mm) 

E  voltage  (V) 

£fev  reversible  voltage  (V) 

Et h  thermoneutral  voltage  (V) 

F  Faraday  constant  (96,485  C  mol- 1 ) 

G  Gibbs  free  energy  (J  N-2  m-3) 

H  rate  of  heat  transfer  (kW) 

J  cell  current  density  (A  m-2) 

/  distance  from  the  manifold  inlet  (m) 

Lc h  flow  channel  length  (m) 

Lm  manifold  length  (m) 

m  mass  flow  rate  (kg  s  - 1 ) 

n  exponent  used  to  represent  frictional/minor  losses 

A  molar  flow  rate  (mol  s  ) 

Aceii  number  of  cells 

Ach  number  of  channels 

N\  number  of  loops  in  a  stack  (Aceii  —  1) 

Ar  rate  of  reactant  consumption  in  the  catalyst  layer 
(mols-1) 

As  number  of  segment  in  a  loop  (6) 

Aturn  number  of  turns  in  flow  channel 
Ao  initially  assumed  molar  flow  rate  (mol  s  ) 

NTU  number  of  heat  transfer  unit 

P  pressure  (Pa) 

Q  volume  flow  rate  (m3  s-1) 

R  universal  gas  constant  (8.3 14  kJ  kmol-1  K-1) 
r  flow  resistance  coefficient 

Re  Reynolds  number 

S  entropy  (JK-1) 

T  temperature  (K) 

7h  temperature  of  heat  source  (K) 

rc, in  inlet  temperature  of  cold  stream  (K) 

U  overall  heat  transfer  coefficient  (W  m-2  K-1) 

V  average  velocity  (ms-1) 

Greek  letters 

a  thermal  diffusivity 

A  difference 

s  effectiveness 

f  stoichiometry 

0  overpotential 

0  flow  direction  convention  (+1  for  clockwise,  —  1 

for  counter  clockwise) 

0  cell  thickness  (m) 


A 

bipolar  plate  effective  height  (m) 

pi 

viscosity  (N  s  m-2) 

P 

density  (kgm-3) 

V 

kinematic  viscosity 

Q 

bipolar  plate  effective  width  (m) 

Subscripts 

A 

anode 

bP 

bipolar  plate 

C 

cathode 

ch 

flow  channel 

cell 

fuel  cell 

drag 

dragged  molecules 

e 

electrode  backing 

f 

friction 

H+ 

protons 

h2 

hydrogen 

i 

loop  number 

in,  inlet 

in/inlet  value 

j 

segment  number 

k 

species  number 

1 

laminar  flow 

loop 

loop 

m 

manifold,  membrane 

max 

maximum  value 

min 

minimum  value 

out,  outlet  out/outlet  value 

02 

oxygen 

ref 

reference 

R 

reacting 

stack 

stack 

t 

turbulent  flow 

trans 

transient  flow  condition 

turn 

total/tum 

w 

water  as  a  coolant 

0 

initial/inlet  value 

tries  (often  for  short  straight  flow  channels)  or  over-simplified 
phenomena  involved.  For  example,  the  entire  PEM  fuel  cell 
(membrane-electrode  assembly)  was  completely  neglected  from 
the  computational  domain  [8],  or  more  often  the  electrochemi¬ 
cal  reactions  in  the  catalyst  layers  leading  to  power  generation  is 
either  neglected  completely  or  assumed  to  have  infinitely  large 
reaction  rate,  hence  the  entire  catalyst  layer,  the  heart  and  soul 
of  fuel  cells,  is  simplified  into  an  infinitesimally  thin  mathemat¬ 
ical  surface — an  unrealistic  and  simplistic  approach.  Because 
of  the  over- simplifications  invoked  in  these  models,  sometimes 
artificial  boundary  conditions  have  to  be  invented,  such  as  the 
fictitious  water  concentration  [7].  Also  due  to  the  significant 
computational  time  required,  these  models  are  often  limited  to 
the  investigation  of  specific  phenomena  occurring  in  the  PEM 
fuel  cell,  such  as  electron  transport  [11,12],  heat  transfer  [13], 
water  transport  [17,18],  etc.  It  is  clear  that  these  CFD-based 
cell  models  are  not  suitable  for  the  analysis  and  modeling  of 
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PEM  fuel  cell  stacks  consisting  of  tens  or  hundreds  of  individual 
cells.  The  third  category  of  the  models  available  in  the  literatures 
[7,24]  take  into  account  both  the  electrochemical  reactions  in 
the  catalyst  layers  (which  are  spatially  resolved  in  the  model), 
and  the  physical  transport  of  reactant  gas  flows,  product  and 
process  water,  heat  and  the  charged  species  in  the  individual 
cells  and  stacks,  hence  it  is  a  fundamentally-based  compre¬ 
hensive  approach,  that  provides  the  predictive  capability  with 
sufficiently  accurate  results.  Further,  the  computational  require¬ 
ment  is  quite  modest,  thus  it  is  selected  for  incorporation  into 
the  present  stack  model. 

In  contrast  with  the  extensive  single  cell  models  available 
in  literature,  quite  limited  number  of  studies  is  available  for 
the  modeling  and  simulation  of  PEM  fuel  cell  stacks.  Thirumu- 
lai  and  White  developed  a  performance  model  of  a  PEM  fuel 
cell  stack  [25]  by  integrating  a  pipe  flow  network  into  the  two- 
dimensional  isothermal  single  cell  model  of  Nguyen  and  White 
[26].  It  has  been  shown  that  the  voltage  variation  is  caused  by 
the  unequal  gas  flow  to  the  individual  cells  in  an  isothermal 
stack.  In  Baschuck  and  Li  [2],  the  pressure  and  mass  flow  rates 
were  obtained  from  a  hydraulic  model  and  was  used  as  operating 
conditions  for  their  cell  model  [27].  The  effect  of  stack  mani¬ 
fold  size  and  the  number  of  gas  flow  channels  per  bipolar  plate 
on  the  variance  of  cell  performances  was  investigated  as  well. 
Karimi  et  al.  [28]  conducted  an  extensive  work  on  optimization 
of  the  stack  configuration  using  flow  network  model  along  with 
the  cell  model  of  Baschuk  and  Li  [27].  The  minor  pressure  loss 
caused  by  complicated  flow  patterns  in  stack  was  shown  to  have 
considerable  impact,  and  the  multiple  inlets  for  the  stack  mani¬ 
fold  were  suggested  to  increase  the  flow  uniformity  of  the  flow 
channels  in  each  cell. 

There  exist  only  few  results  that  deal  with  the  thermal  man¬ 
agement  in  the  fuel  cell  stack.  Maggio  et  al.  [29]  investigated 
the  temperature  and  current  density  distribution  using  a  three- 
dimensional  model.  The  temperature  in  the  cooling  plate,  cool¬ 
ing  water  and  membrane  electrode  assembly  (ME A)  was  found 
through  application  of  conservation  of  energy,  while  the  cell  per¬ 


formance  was  determined  with  the  empirical  relation  of  Patel 
et  al.  [30].  Using  empirical  model  of  Kim  et  al.  [31],  Lee  and 
Lark  carried  out  a  numerical  simulation  for  the  thermo-physical 
behavior  of  fuel  cell  stack  [32].  The  entire  fuel  cell  system 
including  oxidant  supply  and  cooling  pump  was  modeled  in 
Barbir  et  al.  [33]  with  empirical  linear  voltage  relationship. 

In  this  work,  a  non-isothermal  stack  model  has  been  devel¬ 
oped  to  analyze  and  optimize  the  PEM  fuel  cell  stack  perfor¬ 
mance.  The  model  consists  of  a  flow  network  model  for  the 
thermal  and  physical  behavior  of  the  reactant  streams  for  the 
stack  flow  configuration  investigated,  a  thermal  model  for  the 
temperature  distribution  throughout  the  stack  and  a  fuel  cell 
model  based  on  Baschuk  and  Li  [27]  for  the  individual  cell  per¬ 
formance.  The  pressure,  temperature  and  flow  distributions  for 
the  fuel,  oxidant  and  coolant  streams  in  the  stack  are  obtained 
by  incorporating  cell  performance,  heat  generation  and  transport 
phenomena.  The  heat  transfer  to  the  fuel,  oxidant  and  coolant 
streams  are  analyzed  and  the  effects  of  temperature  distribution 
on  the  stack  performance  are  investigated.  The  performance  of 
the  stack  is  evaluated  for  different  stack  design  and  operating 
conditions. 

2.  Mathmatical  models 

2.7.  Stack  flow  model 

Fig.  1  shows  a  schematic  drawing  of  a  PEM  fuel  cell  stack  in 
which  cells  are  combined  in  series  to  generate  sufficient  power 
for  practical  application.  The  individual  cell  consists  of  a  mem¬ 
brane  electrolyte  assembly  (MEA)  sandwiched  between  two 
bipolar  plates.  The  electrode  backings,  catalyst  layers  and  poly¬ 
mer  electrolyte  membrane  are  referred  to  as  MEA  where  direct 
current  is  produced.  In  each  cell,  fuel  and  oxidant  are  supplied 
to  both  sides  of  MEA  through  the  flow  channels  on  the  bipo¬ 
lar  plates.  The  electron  is  generated  in  the  anode  catalyst  layer 
through  oxidation  and  consumed  in  the  cathode  catalyst  layer  by 
reduction.  The  polymer  electrolyte  membrane  conducts  the  pro- 
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®  Coolant  Channel 
©  Bipolar  Plate  (Gas  Flow  Channel) 

©  Electrode  Backing  Layer  -v 

©  Catalyst  Layer 

®  Polymer  Electrolyte  Membrane  Layer 


Membrane  Electrode  Assembly  (MEA) 


Fig.  1.  Schematic  of  a  polymer  electrolyte  membrane  fuel  cell  stack. 
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ton  produced  by  hydrogen  oxidation  to  the  cathode.  Electrode 
backing  and  bipolar  plates  are  conductive  to  deliver  electrons 
from  one  cell  to  the  next.  To  remove  the  heat  produced  by  a 
PEM  fuel  cell,  liquid  water  flows  as  a  coolant  through  cool¬ 
ing  path.  An  additional  plate  is  inserted  to  provide  coolant  path 
between  each  cell  or  the  cooling  path  could  be  integrated  in  the 
bipolar  plate.  Hydrogen  and  oxygen  are  supplied  to  each  cell 
through  the  intake  manifolds  while  excess  reactants  and  reac¬ 
tion  product  exit  to  the  exhaust  manifolds.  The  electrochemical 
reaction  starts  at  the  inlet  of  flow  channel  and  continues  down¬ 
stream  while  the  chemical  energy  is  converted  into  electrical 
energy,  producing  heat  and  water  at  the  same  time.  The  gas  flow 
channels  remove  the  water  within  the  MEA  and  also  supply 
the  humidity  required  to  prevent  polymer  electrolyte  membrane 
dehydration. 

Fig.  2  shows  a  flow  network  for  the  present  model  develop¬ 
ment.  In  order  to  obtain  the  pressure,  stoichiometry,  and  reactant 
composition  for  each  cell,  the  mass  flow  rate  and  pressure  distri¬ 
butions  within  the  stack  must  be  determined;  this  constitutes  the 
stack  flow  model.  Flow  channels  on  bipolar  plate  are  modeled  as 
up  and  downstream  segments  with  one  reaction  site  where  the 
electro-chemical  reaction  and  heat  transfer  occur.  Flow  chan¬ 
nels  are  connected  by  intake  and  exhaust  manifolds  constituting 
a  stack  flow  network.  The  flow  configuration  is  assumed  to  be 
identical  for  the  fuel,  oxidant  and  coolant  streams  in  this  work. 
Only  heat  transfer  occurs  for  the  coolant  stream  while  heat  and 
mass  transfers  occur  simultaneously  for  the  other  two  streams. 
The  pressure,  temperature,  and  compositions  are  assumed  to  be 


uniform  in  each  of  the  loop  segment.  The  inlet  hydrogen  and  oxy¬ 
gen  molar  flow  rates  to  the  anode  and  cathode  are  determined 
according  to  stoichiometry: 


£a  N cell  J Acell 

IF 


^C^cell  *^cell 


4  F 


5 


where  and  fc  are  the  anode  and  cathode  stoichiometries, 
Aceii  the  number  of  cells  in  the  stack,  J  the  current  density,  Acen 
the  active  area  of  the  cell,  and  F  is  the  Faraday  constant.  The 
inlets  of  the  anode  and  cathode  intake  manifolds  are  saturated 
with  water  vapor  for  the  present  simulation.  The  total  flow  rates 
are  obtained  by  adding  the  amount  of  water  vapor  to  that  of 
hydrogen  for  the  anode  stream  and  by  adding  the  amount  of 
water  and  nitrogen  to  that  of  oxygen  for  the  cathode  stream.  The 
mass  must  be  conserved  at  every  node  except  at  the  reaction  site. 
The  rate  of  hydrogen  production  at  the  reaction  site  depends  on 
the  current  density  and  it  is  calculated  using  Faraday’s  law: 


AA^r,h2 


J ^cell 
IF 


where  the  negative  sign  represents  the  fact  that  hydrogen  is 
actually  consumed  at  the  reaction  site.  It  is  assumed  that  the 
exit  stream  of  the  anode  is  fully  saturated  and  the  excess  water 
is  migrated  through  the  polymer  electrolyte  membrane  as  a 
net  result  of  electroosmotic  drag  and  back  diffusion.  The  total 
consumption  rate  in  the  anode  reaction  site  can  be  written  as 


Inlet 


^  L[l][l] 


J[  L[l][4] 
Outlet 

(U-configuration) 


(Intake  Manifolds) 


L[/-l][l]  L[/][l]  L[/+I][l] 


L[AM][1] 


L[AM][4]  j[ 


Outlet 

(Z-configuration) 


///  ire  :  The  amount  of  heat  transfer  at  z'th  loop  (I:  to  the  intake  manifold,  R:  to  the  reaction  site,  E:  to  the  exhaust  manifold) 
L[z][/]  :  A  loop  segment  in  flow  network  (/:  loop  number^’:  loop  segment  number) 

N :  Number  of  cell 

iL 

ANj :  The  amount  of  reaction  at  z  loop  (positive  for  anode,  negative  for  cathode) 

R[i] :  A  reaction  site  at  z',h  loop 


Fig.  2.  A  flow  network  diagram  of  heat  and  mass  transfer  in  a  PEM  fuel  cell  stack;  numbers  (1)— (6)  represent  the  loop  segment  number  at  zth  loop. 
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follows: 


ANa  =  AArh2  +  drag,H20  (4) 

where  Adrag,H2o  is  the  net  water  transport  from  the  anode  to 
the  cathode  through  the  electrolyte  membrane.  In  cathode  cata¬ 
lyst  layer,  oxygen  is  consumed  and  water  molecule  is  produced 
by  the  chemical  reaction.  The  resultant  composition  at  the  exit 
of  reaction  site  as,  considering  the  amount  of  dragged  water 
molecule: 


ANq  =  AAr?q2  +  AA^r?h20  +  Ardrag,H20 


The  production  rate  of  oxygen,  AAr5o2,  and  the  production 
rate  of  water,  A^h2o,  are  given  as 


AA^r,o2 


J Aced 
4  F 


AA^r,h2o 


J  A  QQ 11 

IF 


Water  condensation  may  occur  as  the  amount  of  water 
increases  in  the  cathode.  Hence,  four  species  may  exist  at  the 
cathode  channel  exit,  including  oxygen,  liquid  water,  water 
vapor  and  nitrogen. 

In  Fig.  2,  two  adjacent  gas  flow  channels  and  a  pair  of  intake 
and  exhaust  manifolds  constitute  a  closed  loop  around  a  ME  A. 
The  loop  segment  is  numbered  clockwise  starting  from  the  intake 
manifold  as  shown  in  Fig.  2.  The  pressure  changes  around  the  ith 
closed  loop  should  be  zero  to  satisfy  the  conservation  of  energy, 
i.e. 


Ns 

J2diJAPU  =  0  O'  =  1.  2,  3,  ... ,  N\)  (8) 

j=  i 


where  j  is  a  segment  number,  Ns  the  number  of  segments  in 
a  loop,  Qij  a  sign  convention  for  the  flow  direction,  A Pq  the 
amount  of  pressure  drop  at  the  segment  j  in  ith  loop  and  N\  is  the 
total  number  of  loops  in  a  stack.  In  any  closed  loop,  Oq  is  consid¬ 
ered  to  be  1  for  the  clockwise  and  —  1  for  counterclockwise  flow. 
The  pressure  drop  caused  by  frictional  loss  on  the  inside  wall  of 
flow  channel  is  estimated  based  on  Darcy- Weisbach  equation: 


A  Pi,j  — 


where  Lq  is  the  length  of  loop  segment,  D^jj  the  hydraulic 
diameter,  pij  the  fluid  average  density,  and  Vq  the  flow  speed, 
and  Cfjj  is  the  friction  coefficient  determined  from  the  following 
correlation: 


Ci/ReDhij  Reohi,j  <  2  x  103 
0.316/ Re^jj  4  x  103  <  ReDhiJ  <  105 


(10) 


where  C\  is  a  constant  depending  on  the  flow  path  geometry.  The 
Reynolds  number  in  Eq.  (10)  is  defined  based  on  the  hydraulic 
diameter  of  each  loop  segment  wth  //,/j  representing  average 
viscosity: 


Pi,jVijDh,iJ 
A6 iJ 


The  value  of  Cfjj  is  linearly  interpolated  for  the  value  of 
Reij  between  2000  and  4000.  In  the  anode  and  cathode  streams, 
the  volume  flow  rate  is  obtained  based  on  the  ideal  gas  law  as 
follows: 


Nk'tjRT 

Pk,i,j 


(12) 


where  k  is  the  species  number,  Nkjj  and  Pkjj  are  the  molar 
flow  rate  and  partial  pressure  of  the  species  k,  respectively.  The 
volume  of  the  condensed  liquid  water  is  neglected  in  anode 
and  cathode  downstreams.  The  thermal  properties  of  species 
are  determined  based  on  the  absolute  temperature  and  then  the 
averages  are  calculated  according  to  compositions  following 
Baschuck  and  Li  [2]. 

The  mass  balance  is  applied  to  the  every  node  in  the  intake 
manifold,  reaction  site,  and  exhaust  manifold.  Since  each  closed 
loop  shares  the  vertical  loop  segments  with  adjacent  closed  loop, 
e.g.  L[i\  [2]  is  equivalent  to  L[i  +  1]  [6]  in  Fig.  2,  the  mass  balance 
for  the  vertical  loop  segment  is  expressed  in  terms  of  the  mole 
flow  rate  as 


*i.  t  = 


U  uM— 1,2(+Mnlet) 

°U  1 


(13) 


Nu  5  = 


^•,6^,6  -  A Ni 


6 


i,  5 


NiA  = 


&i-l, 3^-1, 3  -  ^-l,4M-l,4(-^Voutlet) 


0 


i,  4 


(14) 

(15) 


where  A  Ni  is  the  molar  flow  rate  of  reaction,  and  Mnlet  and 
^outlet  are  the  molar  flow  rates  entering  and  leaving  the  fuel 
cell  stack  to  be  considered  only  at  the  stack  inlet  and  outlet 
position.  A  negative  flow  rate  during  iteration  indicates  that 
the  flow  direction  in  a  loop  segment  is  incorrect  and  must  be 
reversed.  The  conservation  of  energy  requires  that  the  summa¬ 
tion  of  pressure  drop  around  any  closed  loop  in  Fig.  2  should 
be  equal  to  zero.  The  pressure  drop  in  a  loop  segment  is  con¬ 
tributed  by  frictional  and  minor  losses,  and  is  related  to  the  flow 
rate: 


A  Pi,  j  —  0i,j[ri,i,jQit  y  T  (.r2,i,j  +  ^3  ,i,j)Qi,j] 


(16) 


where  r\  jj  is  the  coefficient  of  frictional  loss,  the  coeffi¬ 
cient  of  divergence  (confluence)  loss  between  intake  (exhaust) 
manifold  and  flow  channel,  and  is  the  coefficient  of  bending 

loss  in  serpentine  flow  channel.  The  exponent,  uq,  is  dependent 
on  the  flow  condition:  1  for  laminar  flow  and  1.75  for  turbulent 
flow.  r\jj  in  Eq.  (16)  can  be  expressed  for  the  jth  segment  in  ith 
loop,  as  follows: 


rUJ  =  \ 


2-OAfZ)3 


0.158  LijtffpVf 

DlfAl15 

h,i,J  i 


ReDhi,j  <  2  x  103 


4  x  103  <  Reohi,j  <  5  x  105 


(11) 


(17) 
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Correlations  for  V2dj  and  depending  on  the  flow  config¬ 
urations  are  given  in  Table  1  [34] . 

2.2.  Fuel  cell  model 

The  efficiency  of  a  fuel  cell  is  degraded  as  waste  heat  is  pro¬ 
duced  due  to  the  irreversibility  of  electro-chemical  reaction  and 
in  the  transport  of  reactants,  electrons  and  protons  in  the  cell. 
The  amount  of  heat  produced  is  equivalent  to  the  residue  of  the 
maximum  available  chemical  energy  after  producing  electricity. 
In  this  model,  the  amount  of  heat  in  each  cell  is  estimated  from 
the  difference  between  the  thermoneutral  voltage  and  cell  output 
voltage: 

Hi,  total  —  —  is  cell,/)  (18) 

where  Et is  the  thermoneutral  voltage  and  Ece\\j  is  the  cell  volt¬ 
age.  The  thermoneutral  voltage  is  a  property  of  the  fuel  defined 
as 


membrane  are  denoted  by  ijbp,  Oe  and  rjm,  respectively.  The 
voltage  losses  attributed  to  the  bipolar  plate  is  calculated  by 
considering  the  electrode  backing  and  bipolar  plate  as  electrical 
resistances.  The  overpotential  associated  with  the  proton  migra¬ 
tion  in  the  polymer  electrolyte  membrane  is  determined  by  the 
Nernst-Planck  equation  assuming  a  constant  conductivity  for 
the  fully  hydrated  polymer  electrolyte  membrane.  The  anode 
and  cathode  catalyst  layer  overpotentials  are  found  by  con¬ 
sidering  species  conservation,  and  proton  and  electron  migra¬ 
tion  within  the  catalyst  layers.  Proton  and  electron  migration 
through  the  catalyst  layers  are  related  to  the  protonic  and  elec¬ 
trical  current  using  Ohm’s  law.  Species  conservation  requires 
modeling  of  reaction  kinetics  and  mass  transport.  Oxygen  reduc¬ 
tion  is  modeled  with  the  Butler- Volmer  equation  in  the  cath¬ 
ode  catalyst  layer.  The  reversible  voltage,  Erev,  is  calculated 
from  a  modified  Nernst  equation  with  extra  terms  to  account 
for  the  deviation  from  the  standard  reference  temperature  and 
concentration: 


AH 

£rev  =  —  (19) 

nf 

The  individual  cell  voltage  in  the  stack  is  determined  using 
the  single  cell  model  developed  by  Baschuk  and  Li  [27],  and  is 
obtained  from  the  reversible  cell  voltage  Erev  and  the  overpoten¬ 
tials  occurring  in  the  various  components  of  the  PEM  fuel  cell 
as  follows: 

^cell  =  ^rev  On  1 77c  I  ^Obp  2?]e  —  T]m  (20) 

where  rja  and  rjc  are  the  overpotentials  attributed  to  the  anode  and 
cathode  catalyst  layers,  respectively.  The  voltage  losses  caused 
by  the  bipolar  plate,  electrode  backing  and  polymer  electrolyte 


(21) 


where  AG  is  the  change  in  Gibbs  free  energy,  AS  the  change  in 
entropy,  R  the  universal  gas  constant,  T  the  absolute  temperature 
of  the  MEA,  and  Q  denotes  the  concentration  of  species  i  at  the 
electrode  backing/catalyst  layer  in  mol  m-3 .  The  reference  value 
of  temperature,  Tref,  is  298.15  K  while  the  reference  concentra¬ 
tions  for  hydrogen  and  oxygen  are  22.22  and  7.033  mol  m-3 
respectively. The  reference  concentrations  represent  the  concen- 


Table  1 

Branch,  confluence,  and  bend  head  loss  correlations  [34] 


Configuration 

Diverging  T 


i 


J 

1 


3 

Converging  T 


Bending 


l 


Correlations 


A  Pi 


j 


n,  V2 

=  k i-j  fi~Y  (j  =  2,  3)  ki-2 


a 1  +a2  (^)  +ai(J±)1 


2 


APhj  =  k1-J^f(j  =  2,  3)  ki-2 


al  +02  (Jf)  +a3(^)2 


2 


2 

APbend  =  k^L-  k  =  k\oc  +  kfr  k\oc  =  — °'21025  (90°  bend,  square  duct)  kfr  =  0.0175  ^  SX  (S  in  degrees)  k  = 

z  (7?0/«h)  '  “h 

/  / - \  -2/3 

32  \^Redh  y  ^  J  (laminar) 
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tration  of  hydrogen  and  oxygen  in  Nafion  at  a  temperature  of 
298.15  K  and  partial  pressure  of  1  atm. 

2.3.  Heat  transfer  model 

Within  the  framework  of  the  present  analysis,  each  cell  can 
be  considered  as  a  parallel  flow  heat  exchanger  in  which  heat 
is  transferred  from  the  MEA  to  the  fuel,  oxidant  and  coolant 
streams.  The  pressure,  temperature  and  flow  rates  of  the  streams 
must  be  obtained  incorporating  individual  cell  operation  since 
the  thermal  and  physical  flow  conditions  and  the  individual  cell 
performance  are  dependent  on  each  other.  In  this  work,  the 
effectiveness-NTU  method  is  applied  to  determine  the  amount 
of  heat  transfer  and  the  temperature  of  the  each  loop  segment. 
The  log  mean  temperature  difference  (LMTD)  method  is  not 
eligible  for  the  present  study  since  the  unknown  values  are  not 
only  the  outlet  temperature  but  also  the  temperature  of  the  heat 
source  (i.e.  the  temperature  of  MEA).  The  analysis  is  based  on 
the  following  assumptions: 

1.  The  individual  MEA  has  a  uniform  temperature. 

2.  The  reactant  flow  and  thermal  conditions  are  fully  developed, 
neglecting  the  entrance  effect. 

3.  Thermal  resistances  at  the  surface  of  the  flow  channel  includ¬ 
ing  fouling  effects  are  negligible. 

4.  The  convective  heat  transfer  on  the  outside  surface  of  the  fuel 
cell  stack  is  negligible. 

The  number  of  heat  transfer  unit  (NTU)  is  a  dimensionless 
parameter  that  is  widely  used  for  the  heat  exchanger  analysis 


and  is  defined  as 
UA 

NTU  =  -  (22) 

Cmin 

where  U  is  the  overall  heat  transfer  coefficient,  A  the  area  of 
heat  transfer  and  Cm jn  is  the  minimum  heat  capacity  rate.  The 
overall  heat  transfer  coefficient  is  equivalent  to  the  heat  transfer 
coefficient  h  since  conduction  and  other  thermal  resistances  are 
neglected.  The  value  of  h  is  obtained  from 


(23) 


and  the  Nusselt  number  is  determined  from  the  following  cor¬ 
relations  for  the  internal  flow  [35]: 


Nuoh  = 


3.66  ReDh  <  2000,  Pr  >  0.6 

0.023  Re4f  Pr0A  ReDh  >  10000,  Pr  >  0.6) 

(24) 


where  Pr  is  the  Prandtl  number  in  each  loop  segment,  defined 
as: 

Pr=-  (25) 

a 

where  v  is  the  kinematic  viscosity  and  a  is  the  thermal  diffusiv- 
ity.  The  Nusselt  number  of  3.66  for  the  laminar  flow  is  based 
on  constant  surface  temperature,  consistent  with  the  assumption 
of  the  uniform  temperature  in  MEA  [35].  The  Nusselt  number 
within  the  transition  flow  condition,  2000  <  Reph  <  10,  000,  is 
linearly  interpolated.  The  value  of  the  effectiveness,  s,  in  each 


(a)  U-configuration 


(b)  Z-configuration 


Fig.  3.  Schematic  of  temperature  distribution  and  heat  transfer  in  the  flow  network. 
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flow  loop  segment  can  be  found  according  to  the  following  cor¬ 
relation  against  the  constant  temperature  heat  source: 

s  —  l  —  exp(-NTU)  (26) 


following  equations  are  obtained  for  the  U  configuration  stack: 


7/,R,in  —  T'i,  in 

7/,/, in  —  Ti  in 


(28) 

(29) 


The  amount  of  heat  transfer  is  proportional  to  the  effec¬ 
tiveness,  the  heat  capacity  rate,  and  the  temperature  difference 
between  the  ME  A  and  any  other  flow: 

H  =  ^Cmin(7h  —  Tc{n )  (27) 

where  Cm jn  is  the  minimum  heat  capacity  rate,  7h  the  tempera¬ 
ture  of  the  heat  source  and  Tc^n  is  the  inlet  fluid  temperature  of 
the  cold  stream,  respectively. 

Fig.  3  is  given  to  elucidate  the  mathematical  scheme  for  the 
temperature  distribution  and  the  heat  transfer  in  U  and  Z  con¬ 
figuration.  The  temperature  of  the  ME  A  and  the  amount  of  heat 
transfer  to  the  intake  manifold,  reaction  site,  and  exhaust  man¬ 
ifold  is  determined  from  the  given  inlet  temperature  and  the 
total  amount  of  heat  generated,  //total,/-  At  steady  sate,  the  heat 
generated  in  each  MEA  is  discharged  completely  so  that  the 
temperature  of  the  MEA  will  remain  unchanged.  By  applying 
the  energy  and  mass  conservation  at  each  node  in  Fig.  3(a),  the 


7/,R,out 


Z/,R,in  T  . 


H 


i,  R 


Pli,R^  p,i,R 


(30) 


7/,  out  — 


^/,R^/?,/,R7/,R,out  H-  ^/,E,out^/?,/,E7/,E,out 
Pli,  out  C  p,i,  out 


(31) 


The  amount  of  heat  transfer  rate  to  the  intake  manifold,  reac¬ 
tion  site,  and  exhaust  manifold  can  be  estimated  as  follows: 


Hi, I  =  Si,imi,iCp,i,i(Ti,MEA  ~  //, I, in)  (32) 

Hi, R  =  £/,R^/,R<7 p,/,R(7i‘,MEA  -  7},R?in)  (33) 

77/,e  =  £/,e^/,eC /?,/,e(7/,mea  -  Ti,E,in)  (34) 

The  total  amount  of  heat  should  be  equivalent  to  the  summa¬ 
tion  of  above  three  terms  as 


Hi,  i  +  77/, r  +  /7/,e  =  Hi,  total  (35) 

The  unknown  variables  are  TiXin,  7),R,m,  TiXou t,  Hx j,  HUR, 
7//,e,  7^/, out,  7/>e> out  and  7),  mea  inEqs.  (28)— (35).  The  temperature 


START 


Qa 


T Temperature 
P:  Pressure 

C:  Species  concentration 
Q :  Volume  flow  rate 
E:  Voltage 

Ev:  Convergence  criterion  of  voltage 
Q.  Stoichiometry 


Qc 


Tceib  Pf  Qa,  Qc,  Ca,  Cc 


END 


Fig.  4.  A  flow  chart  of  the  stack  model. 


452 


J.  Park,  X.  Li  /Journal  of  Power  Sources  162  (2006)  444-459 


distribution  and  the  amount  of  heat  transfer  are  determined  from 
the  first  cell.  For  the  first  cell  in  U  configuration  stack,  7)?0ut 
corresponds  to  the  temperature  of  the  coolant  out  of  the  stack, 
therefore  it  can  be  estimated  from  the  total  heat  generation  of 
the  stack  as  follows: 


E^total 

IP  total  C  p  ,out 


^l,in  + 


(36) 


Similarly  for  Z  configuration  stack,  Eq.  (31)  is  replaced  by 


rj,  ™i,RC  i  ,E,out  Cpj-l^Ti-l  ,E,0ut  /q-7\ 

1  i, E,in  = - ; - — -  (3  /) 

m,E^p,i,E 

where  7)?E,in  equals  to  7);R;0ut  for  the  first  cell  since  7)_i?out  does 
not  exist  at  the  first  cell  for  Z  configuration.  A  set  of  linear 
Eqs.  (28)-(35)  can  be  solved  analytically  or  an  iterative  solution 
can  be  found  by  modifying  the  temperature  of  the  MEA  until 
Eq.  (35)  is  satisfied.  The  mathematical  scheme  for  the  anode 
and  cathode  gas  streams  is  similar  to  that  of  the  coolant  stream. 
For  fast  convergence,  the  heat  transfer  to  the  anode  and  cathode 
streams  is  neglected  for  the  first  iteration.  Since  coolant  water  has 
much  larger  heat  capacity  compared  to  the  other  two-reactant  gas 
streams,  most  of  heat  generated  in  the  MEA  should  be  transferred 
to  the  coolant  stream.  The  amount  of  heat  transfer  to  the  anode 


Table  2 


Parameters  and  properties  used  in  the  present  PEM  fuel  cell  stack  simulations 


Component 

Parameter 

Value 

Bipolar  plate 

Width  (Q) 

12  x  10-2  m 

Height  (A) 

12  x  10-2  m 

Thickness  (®) 

0.012m 

Gas  flow  channel 

Channel  length  (/ch) 

3.0  x  10“ 1  m 

(anode,  cathode) 

Number  of  channel  (Ach) 

10 

Number  of  U-turn  ( Nt ) 

3-8 

Hydraulic  diameter  (df) 

0.75-2.25  mm 

Manifold  hydraulic  diameter 

(AO 

12.5-20.0  mm 

Gas  inlet  temperature  ( Tm ) 

353. 15K 

Coolant  path 

Path  length  (/ch,w) 

3.0  x  10-1  m 

Hydraulic  diameter  (dh,w) 

1.5  mm 

Number  of  channel  (Nch,w) 

10 

Number  of  turn  (/VUw) 

3 

Manifold  hydraulic  diameter 
(Ai,w) 

20.0  mm 

Coolant  inlet  temperature 

(7in,w) 

353. 15K 

Coolant  flow  rate  (<2W) 

0.85  x  10"4  to 

2.84  x  10-3  m3  s-1 

Reynolds  number  (Re of) 

1.48  x  104  to 

4.95  x  105 

Stack 

Number  of  cell  (Vce n) 

51 

Anode  stoichiometry  (£a) 

1.2 

Cathode  stoichiometry  (£c) 

2.0 

Current  density  (J) 

0.5  As-2 

Electrode  backing 

Thickness  (5) 

2.5  x  10-4  m 

Porosity  (0) 

0.4 

Catalyst  layer 

Thickness  ( 8C ) 

2.0465  x  10"5  m 

Membrane 

Thickness  (<5m) 

1.64  x  10-4  m 

and  cathode  is  estimated  by  Eqs.  (32)-(34)  considering  the  heat 
of  water  condensation. 


3.  Numerical  procedures 


Fig.  4  shows  a  flow  diagram  of  the  present  stack  model.  The 
stack  model  consists  of  three  major  modules:  a  flow  solver,  a  heat 
transfer  solver,  and  a  fuel  cell  module.  In  the  flow  solver,  pressure 
and  flow  rates  of  species  for  the  anode,  cathode,  and  coolant 
streams  are  obtained  according  to  the  given  inlet  conditions.  The 
flow  solver  is  based  on  the  Hardy  Cross  method  which  needs  an 
iterative  procedure  to  find  pressure  and  flow  rate  distributions 
in  a  flow  network  [36].  Starting  from  an  assumed  flow  rate  in 
each  segment  of  all  the  loops,  the  flow  rate  of  each  loop  segment 
is  modified  observing  the  mass  and  energy  conservation  until 
a  converged  flow  distribution  is  found.  With  Q/j,prev  being  the 
volume  flow  rate  of  the  previous  iterative  step,  the  flow  rate  is 
corrected  as 


Qij  —  Q  ij,  prev  +  A  Qi 


(38) 


Fig.  5.  Typical  temperature  distributions  of  the  MEA  and  coolant  streams 
along  the  intake  and  exhaust  manifold  as  well  as  the  average  tempera¬ 
ture  in  the  flow  channel  for  a  fuel  cell  stack;  D^  =  20mm,  <4  =  0.75  mm, 
<2W  =  1.42  x  10-4  m3  s-1 .  (a)  U  configuration  and  (b)  Z  configuration. 
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Substituting  into  Eq.  (16)  for  Qtj  and  keeping  the  terms  up 
to  the  first  order  about  A Qi  yields: 


~  @i,jrl,i,jQ i, j2prev  +  @ i,j(r2,i,j  +  ^3,i',j)Gi,j,prev 
AQi[nijr\jj  +  (?2  ,i,j  +  r3  ,i,j)] 


(39) 


Substituting  this  expression  into  Eq.  (8)  we  have 

XI  j=l@i,j(rl,i,j  2/j(prev+  (r2 ,i,j+  r3,i,j)Qi,j,prev) 


AQi  = 


^j=l^ij(niJrl,ijQJijjw>rly  +  (r2JJ  +  r3,  ij)) 


(40) 


Then  Eqs.  (13)— (17),  (38)  and  (40)  are  solved  iteratively  until 
the  summation  of  Eq.  (39)  around  a  closed  loop  i  converges 
within  a  preset  tolerance,  which  is  less  than  10-3  Pa  in  every 
loop  segment. 

After  reaching  a  converged  solution  for  flow  solver,  the 
amount  of  heat  transfer  to  each  stream  is  determined  in  asso¬ 
ciation  with  the  flow  rates  and  heat  generation  in  each  cell  as 
described  in  the  previous  section.  Then  temperature,  pressure 
and  species  concentrations  are  sent  to  the  fuel  cell  model  to 
obtain  output  voltage  of  each  cell.  All  material  properties  are 


Fig.  6.  Temperature  distributions  of  the  gas  streams  in  anode  and  cathode  stack 
manifolds;  D h  =  20  mm,  d h  =  0.75  mm,  Qw  =  1 .42  x  10~4  m3  s_  1 .  (a)  U  config¬ 
uration  and  (b)  Z  configuration. 


renewed  according  to  the  temperature  and  pressure  for  the  next 
iteration.  The  iteration  is  carried  out  until  the  maximum  differ¬ 
ence  of  the  each  cell  voltage  from  the  value  of  previous  iteration 
is  less  than  10-5  V. 

4.  Results  and  discussion 

In  the  present  study,  analysis  has  been  carried  out  for  a  PEM 
fuel  cell  stack  which  consists  of  51  cells  with  144  cm2  active 
cell  area.  The  current  density  is  fixed  at  0.5  A  cm-2  generating 
about  2  kW  as  output  power.  Table  2  lists  the  dimensions,  prop¬ 
erties,  and  operating  conditions  in  the  present  study.  A  total  of 
five  and  seven  different  sizes  for  the  manifold  and  flow  channel 
designs  with  nine  different  coolant  flow  rates  have  been  con¬ 
sidered.  This  resulted  into  a  total  of  630  cases  investigated,  or 
315  cases  for  each  of  the  U  and  Z  stack  configurations.  For 
each  case,  it  took  about  30  min  of  computation  time  in  a  Linux 
PC  cluster  built  by  10  P-4  CPUs  with  20  GB  of  main  mem¬ 
ory.  This  indicates  that  the  present  numerical  algorithm  is  very 


Fig.  7.  Voltage  distributions  in  U  configuration  fuel  cell  stack  for  vari¬ 
ous  coolant  flow  rates:  <2w  =  0.85  x  10_4-2.84  x  10-3  m3  s_1  (Reoh  =  1.48  x 
104-4.95  x  105).  (a)  Design  A:  Dh  =  20mm  and  dh  =  0-75mm  and  (b)  design 
B:  Dh  =  12.5  mm  and  d h  =  2.0  mm. 


454 


J.  Park,  X.  Li  /Journal  of  Power  Sources  162  (2006)  444-459 


computationally  efficient  for  a  stack  of  this  size,  and  is  useful 
for  practical  stack  design  calculations.  Section  4.1  presents  the 
flow  and  temperature  distributions  for  a  typical  operating  con¬ 
dition  and  their  effects  on  stack  performance.  Then,  in  Section 
4.2,  stack  optimization  study  is  presented  based  on  the  net  stack 
power. 

4.1.  Effect  of  flow  and  temperature  distribution 

Flow  non-uniformity  among  the  cells  in  a  stack  is  known  to 
induce  performance  degradation  of  fuel  cell  stacks  [26,28].  Suf¬ 
ficient  flow  uniformity  can  be  ensured  by  enlarging  manifold 
diameter  so  that  the  manifold  behaves  like  a  plenum,  however  it 
is  limited  by  the  size  of  bipolar  plate.  The  individual  cell  perfor¬ 
mance  also  depends  on  the  cell  temperature  which  is  controlled 
by  the  coolant  flow  rate  in  a  fuel  cell  stack. 

Analysis  has  been  carried  out  for  two  common  stack  designs: 
U  and  Z  configurations.  The  inlet  and  outlet  of  the  stack  manifold 
are  located  at  the  first  cell  in  case  of  the  U  configuration  while 
the  outlet  of  Z  configuration  is  located  at  the  last  cell.  A  typ¬ 
ical  temperature  distribution  of  the  MEA  and  coolant  streams 
for  the  U  and  Z  configuration  stacks  is  shown  in  Fig.  5.  The 
inlet  temperature  of  the  coolant,  anode  and  cathode  streams  are 
fixed  at  353.15  K  (80  °C)  in  this  study.  The  temperature  of  the 
coolant  stream  gradually  increases  along  the  intake  manifold. 
The  majority  of  the  temperature  rise  occurs  as  the  coolant  passes 
through  the  flow  channels  on  the  cooling  plate,  due  to  the  heat 
transfer  from  the  MEA.  In  the  exhaust  manifold,  the  tempera¬ 
ture  for  the  U  stack  configuration  seems  to  decrease  along  the 
flow  direction  since  the  coolant  temperature  out  of  each  cell  is 
higher  at  upstream  of  the  exhaust  manifold.  The  temperature 
of  the  coolant  stream  can  exceed  that  of  the  MEA  due  to  the 
heat  transfer  and  increasing  mass  flow  rate  along  the  exhaust 
manifold.  The  temperature  distribution  for  the  Z  configuration, 
however,  increases  along  the  flow  direction,  with  almost  the 
same  outlet  temperature.  A  comparison  reveals  that  the  MEA 
temperature  for  the  U  configuration  is  slightly  higher  than  that 
for  the  Z  configuration.  Because  of  the  difference  for  the  flow 
direction  in  the  exhaust  manifold,  U  configuration  results  in  the 
temperature  in  the  exhaust  manifold  being  higher  than  the  MEA 
for  the  majority  of  the  cells  in  the  stack;  while  for  the  Z  con¬ 
figuration  the  temperature  in  the  exhaust  is  always  less  than  the 
MEA  temperature. 

The  temperature  distributions  of  the  gas  streams  in  the  anode 
and  cathode  stack  manifold  are  similar  to  those  of  the  coolant 
stream,  as  shown  in  Fig.  6.  The  rate  of  the  temperature  increment 
of  the  cathode  stream  is  higher  for  a  significant  number  of  cells  in 
the  stack  than  that  of  the  anode  stream  since  the  cathode  stream 
is  turbulent  at  the  inlet  of  the  intake  manifold  resulting  in  a  larger 
Nusselt  number.  However,  near  the  end  of  the  intake  manifold, 
the  temperature  of  the  anode  gas  stream  increases  faster  com¬ 
pared  to  the  cathode  stream,  because  of  rapidly  decreasing  mass 
flow  rate  of  the  anode  stream  there.  The  amount  of  heat  transfer 
to  the  gas  streams  is  found  to  be  3.2  W  for  the  anode  stream  and 
20.4  W  for  the  cathode  stream  while  the  total  amount  of  heat 
generation  is  3480.5  W,  compared  to  the  stack  electrical  power 
output  of  1916.1  W. 


Fig.  7  shows  the  cell  voltage  distributions  in  two  distinct 
designs  of  U  configuration  as  a  function  of  the  coolant  flow 
rates.  In  Fig.  7(a),  the  cell  voltage  increases  along  the  stack 
while  it  decreases  in  Fig.  7(b).  This  trend  is  more  apparent 
for  lower  coolant  flow  rates.  For  the  stack  design  A,  it  has  a 
small  flow  channel  diameter  (<4  =  0.75  mm)  and  relatively  large 
manifold  diameter  ( D h  =  20  mm),  and  is  hence  characterized  by 
high  flow  uniformity.  On  the  other  hand,  the  flow  uniformity 
is  low  for  the  stack  design  B  since  it  has  a  large  flow  chan¬ 
nel  diameter  ( d h  =  2  mm)  and  relatively  small  manifold  diameter 
(Dh  =  12.5  mm).  The  temperature  effect  is  dominant  in  Fig.  7(a) 
due  to  high  flow  uniformity  and  this  trend  is  similar  for  the  Z 
configuration  shown  in  Fig.  8(a).  The  flow  transition  from  turbu¬ 
lent  to  laminar  flow  is  responsible  for  the  change  in  the  slope  of 
the  cell  voltage  and  temperature  profiles  for  cell  numbers  above 
40  shown  in  Fig.  7(a).  The  voltage  decrement  in  Fig.  7(b)  is 
attributed  to  the  effect  of  flow  variances  along  the  intake  mani¬ 
fold  caused  by  the  lower  flow  uniformity  of  design  B. 

Different  from  the  design  B  for  the  U  configuration,  the  volt¬ 
age  profile  of  the  design  B  for  the  Z  configuration  shows  a  local 
minimum  around  the  middle  of  the  stack  by  the  combined  effect 
of  temperature  and  concentrations.  This  phenomenon  is  caused 
by  the  different  pressure  distribution  for  the  Z  configuration  as 


Fig.  8.  Voltage  distributions  in  Z  configuration  fuel  cell  stack  for  various  coolant 
flow  rates:  (a)  D h  =  20  mm,  d h  =  0.75  mm,  (b)  =  12.5  mm,  d h  =  2.0  mm. 
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Fig.  9.  Pressure  distributions  in  anode  and  cathode  manifolds;  (a)  anode  (U  configuration),  (b)  anode  (Z  configuration),  (c)  cathode  (U  configuration),  (d)  cathode 
(Z  configuration);  Qw  =  1.42  x  10~4  m3  s_1;  design  A:  =  20  mm  and  <4  =  0.75  mm;  design  B:  =  12.5  mm  and  =  2.0  mm. 


shown  in  Fig.  9.  High  flow  uniformity  is  characterized  by  high 
average  pressure  and  low  pressure  drop  in  the  intake  manifold, 
resulting  in  low  concentration  variation  in  anode  and  cathode 
flow  channels  for  design  A.  However  a  different  distribution  is 
observed  in  the  exhaust  manifold  due  to  different  exit  positions, 
i.e.  the  pressure  increases  with  the  cell  number  for  the  U  con¬ 
figuration  while  decreases  for  the  Z  configuration.  This  results 
in  the  pressure  difference  having  a  minimum  near  the  middle 
of  the  Z  configuration,  which  is  more  evident  for  the  cathode 
manifolds. 

To  clarify  the  effect  of  the  flow  and  temperature  distribu¬ 
tions  on  the  stack  performance,  the  distributions  of  temperature, 
species  concentration,  and  stoichiometry  in  the  anode  and  cath¬ 
ode  manifolds  are  shown  in  Fig.  10.  The  MEA  temperature  is 
also  affected  by  the  uniformity  of  the  flow  in  the  gas  flow  chan¬ 
nels.  The  small  differences  in  the  cell  inlet  temperature  T^n 
caused  by  the  degradation  of  flow  rates  in  the  anode  and  cath¬ 
ode  flow  channels  are  accumulated  along  the  intake  manifold, 
resulting  into  a  clear  temperature  difference  for  the  last  few  cells. 
The  pressure  distribution  in  the  anode  and  cathode  streams  are 
considerable  for  the  design  B  due  to  the  lower  flow  uniformity, 
which  leads  to  a  considerable  concentration  variations  for  design 
B  in  Fig.  10.  As  a  result,  the  cell  voltage  distributions  for  design 
A  in  Figs.  7  and  8  resemble  the  temperature  distributions  in 
Fig.  10  when  the  effects  of  concentration  distribution  are  dimin¬ 
ished.  On  the  other  hand,  voltage  distributions  for  design  B  in 
Figs.  7  and  8  are  similar  to  that  of  oxygen  concentration  in  the 
cathode  in  Fig.  10.  The  effect  of  concentration  degradation  in 


cathode  is  dominant  since  the  majority  of  the  voltage  loss  occurs 
in  cathode  [27] . 

From  the  above  discussion,  it  becomes  clear  that  a  uni¬ 
form  voltage  distribution  for  the  cells  in  a  stack  could  be 
achieved  without  a  uniform  distribution  of  the  reactant  flows 
and  temperature — a  condition  difficult  to  achieve.  For  example, 
by  examining  the  effect  of  flow  and  temperature  distribution, 
as  shown  in  Fig.  7  for  the  U  stack  configuration,  the  tempera¬ 
ture  distribution  will  increase  the  cell  voltage  along  the  stack, 
while  the  flow  non-uniformity  will  cause  the  cell  voltage  to 
decrease  along  the  stack.  Then  a  judicially  selected  distribution 
of  the  reactant  flows  and  temperature  along  the  stack  could  bal¬ 
ance  their  respective  effect  so  that  a  uniform  cell  voltage  might 
be  achieved.  This  provides  a  new  approach  to  optimal  stack 
design. 

4.2.  Optimization  of  manifolds  and  flow  channel  designs 

In  the  previous  section  it  is  shown  that  the  cell  voltage  is 
affected  by  the  distribution  of  the  reactant  gas  flows  and  temper¬ 
ature.  The  effect  of  flow  variation  can  be  minimized  by  ensuring 
sufficient  flow  uniformity  through  an  enlarged  manifold  and 
reduced  flow  channel  sizes.  However  excessively  small  flow 
channels  can  lead  to  excessive  pumping  power  required  to  drive 
the  reactant  gas  flow,  in  addition  to  other  problems  such  as  manu¬ 
facturing  difficulty.  The  cell  voltage  is  increased  with  increasing 
temperature  but  a  PEM  fuel  cell  stack  cannot  be  operated  at  much 
higher  than  80  °C.  Optimal  performance  of  a  PEM  fuel  cell  stack 
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Fig.  10.  ME  A  temperature,  species  concentration  and  stoichiometry  distributions  in  anode  and  cathode  manifolds;  (a)  U  configuration  and  (b)  Z  configuration; 
<2W  =  1.42  x  10-4  m3  s_1;  design  A:  Dh  =  20  mm  and  d h  =  0.75  mm;  design  B:  Dh  =  12.5  mm  and  dh  =  2.0  mm. 


can  be  achieved  by  optimizing  the  size  of  the  manifold  and  flow 
channel  as  well  as  the  coolant  flow  rate. 

Fig.  1 1  presents  the  effect  of  coolant  flow  rate  on  the  aver¬ 
age  MEA  temperature,  average  cell  voltage,  pumping  power  and 
net  stack  output  power  for  the  designs  A  and  B.  The  pumping 
power  is  calculated  as  the  volumetric  flow  rate  times  the  amount 
of  pressure  drop  in  each  stream.  The  average  temperature  of 
the  MEA  decreases  as  the  coolant  flow  rate  increases  and  it  is 
likely  to  converge  to  the  inlet  coolant  temperature.  The  effect 
of  stack  design  on  the  average  temperature  of  MEA  is  observed 
although  it  becomes  smaller  with  increased  coolant  flow  rate. 
The  high  flow  uniformity  with  higher  species  concentrations  in 
design  A  results  in  a  higher  cell  voltage  compared  to  that  of 
design  B  for  both  U  and  Z  configurations.  However  the  design 
A  demands  more  pumping  power  especially  for  the  cathode  gas 
streams.  The  total  net  pumping  power  can  be  found  accord¬ 
ing  to  the  coolant  flow  rate  since  the  pumping  powers  for  the 


anode  and  cathode  are  already  determined  from  the  inlet  stoi¬ 
chiometries  and  total  pressure  drop  for  the  given  configuration. 
The  reasonable  range  for  the  coolant  flow  rate  may  vary  from 
0.85  x  10  4  to  3.0  x  10  3  m3  s  1 .  Maximum  temperature  of  the 
MEA  reaches  373.15  K  (=100  °C)  when  the  coolant  flow  rate  is 
less  than  0.8  x  10  4m3s  1 . 

Fig.  12  shows  a  three  dimensional  plot  of  the  net  stack  power 
for  wide  ranges  of  manifold  and  flow  channel  diameters.  The 
individual  surface  plot  represents  the  result  for  a  given  coolant 
flow  rate  and  it  exhibits  major  dependency  on  the  flow  chan¬ 
nel  diameter.  The  variance  of  the  net  power  with  respect  to  the 
flow  channel  diameter  becomes  less  considerable  for  the  larger 
manifold  diameter.  For  each  cooling  flow  rate,  there  exists  an 
optimal  flow  channel  diameter  for  which  the  net  stack  power 
achieves  a  maximum  value.  Low  flow  uniformity  is  responsi¬ 
ble  for  the  power  degradation  at  large  flow  channel  diameters 
while,  at  small  flow  channel  diameters,  the  net  stack  power  is 
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Fig.  11.  Average  MEA  temperature,  average  cell  voltage,  total  pumping  power  and  net  stack  power  according  to  various  coolant  flow  rate,  (a)  U  configuration  and 
(b)  Z  configuration;  Qw  =  1.42  x  10“4  m3  s-1 ;  design  A:  D h  =  20  mm  and  d h  =  0.75  mm,  design  B:  =  12.5  mm  and  d h  =  2.0  mm. 


deteriorated  by  the  excessive  pumping  power.  It  is  also  seen  that 
the  stack  power  increases  as  the  coolant  flow  rate  is  decreased 
because  of  higher  average  cell  temperature.  In  Fig.  13,  the  net 
stack  power  is  compared  between  the  U  and  Z  configurations 


for  various  coolant  flow  rates.  The  maximum  net  stack  power  is 
obtained  when  the  manifold  and  flow  channel  diameters  are  20 
and  1.5  mm  (design  A  in  Fig.  13),  respectively.  Although  it  is 
not  considerable,  the  maximum  net  power  of  U  configuration  is 
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Fig .  1 2 .  A  three-dimensional  plot  of  the  net  stack  power  for  various  stack  designs 
and  coolant  flow  rates  (Reoh  =  1-48  x  104-4.95  x  105).Z4  =  12.5-20  mm,  and 
<4  =  0.75-2.25  mm.  (a)  U  configuration  and  (b)  Z  configuration. 


Coolant  Flow  Rate  (xlO'WVs) 

Fig.  13.  Comparison  of  the  net  stack  power  between  U  and  Z  configurations  for 
various  coolant  flow  rates;  design  A:  Z4  =  20  mm  and  <4  =  1.5  mm;  design  B: 
Dh  =  1 .25  mm  and  <4  =  2.0  mm. 


slightly  higher  than  that  of  Z  configuration  for  the  entire  range 
of  the  coolant  flow  rates  investigated.  On  the  other  hand,  the  net 
power  for  the  U  configuration  is  reduced  considerably  when  the 
flow  uniformity  is  poor  such  as  for  design  B  shown  in  Fig.  13. 
The  effect  of  flow  variance  on  net  stack  power  is  much  less  for 
the  Z  configuration  since  it  has  higher  flow  uniformity  than  the  U 
configuration  for  the  same  inlet  conditions  as  shown  in  Fig.  10. 
The  Z  configuration  is  hence  preferable  when  the  sufficient  flow 
uniformity  cannot  be  ensured. 

5.  Conclusions 

In  this  study,  a  non-isothermal  stack  model  has  been  devel¬ 
oped  for  PEM  fuel  cells,  and  the  effect  of  reactant  flows  and 
temperature  distribution  on  the  stack  performance  has  been  ana¬ 
lyzed.  It  is  shown  that  the  flow  and  temperature  distribution  have 
a  different  influence,  and  a  judicial  matching  of  their  distribution 
can  provide  the  ideal  uniform  cell  voltage  distribution.  Optimal 
stack  design  for  the  maximum  net  stack  power  has  also  been 
investigated  through  minimizing  performance  degradation  and 
excessive  pumping  power.  The  effect  of  flow  distribution  on  the 
stack  performance  is  found  to  be  considerably  less  for  the  Z 
configuration. 
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